Bold Line Diagrammatic Monte Carlo Method: General formulation and application 
to expansion around the Non-Crossing Approximation 



Emanuel GuU,^ David R. Reichman,^ and Andrew J. Millis^ 

^Department of Physics, Columbia University, New York, New York 10027, USA 
'^Department of Chemistry, Columbia University, New York, New York 10027, USA 

(Dated: October 4, 2010) 

We present a general framework for performing "bold-line" diagrammatic Monte Carlo calcula- 
tions using an analytical partial resummation as a starting point for a stochastic summation of all 
diagrams. As a stringent test case we assess the accuracy of the method by solving the equations of 
single-site dynamical mean-field theory, using the non-crossing approximation as a starting point. 
We establish the validity of the starting approximations and show that the bold method provides a 
very accurate treatment of the Mott-insulating phase. 



I. INTRODUCTION 

Observables in quantum field theories may be ex- 
pressed as infinite sums of Feynman diagrams. During 
the last decade it has been realized that stochastic meth- 
ods may be used to estimate the entire diagrammatic 
seriesiii^ This requires a diagram-generating procedure 
which is ergodic (all diagrams contributing to the series 
must be generated) and which is such that each diagram 
is generated with a probability proportional to its weight 
in the series to be studied. These problems have been 
solved^Ti^ and stochastic diagrammatic methods currently 
provide the best estimates for properties of the moder- 
ately correlated regime of the two-dimensional Hubbard 
model as well as successful impurity solvers for dynamical 
mean-field (DMFT)^ and non-equilibrium^i^ problems. 

Analytical studies over many years have established 
that while direct evaluation of low-order terms in a 
perturbation theory is rarely reliable, partial resumma- 
tions of infinite series of diagrams often capture much 
of the relevant physics. Partial resummation replaces 
bare propagators (typically denoted by light lines in dia- 
grams) with rcnormalized propagators (typically denoted 
by heavy or bold-face lines); thus, the diagrams involv- 
ing further corrections to a partial resummation are often 
referred to as bold-line diagrams. Given these successes, 
it is natural to ask if stochastic techniques can be used 
to estimate the further corrections to an infinite partial 
resummation or, in other words, to develop a bold-line 
numerics. References lol and [lOl present one such method, 
a Monte Carlo computation based on the expansion of 
the self-energy of the polaron problem (single particle 
coupled to an oscillator bath) about the ladder resumma- 
tion; however, a bold-line method for a truly many-body 
problem has heretofore not existed. 

In this paper we present a method for the stochas- 
tic evaluation of a many-body fermionic bold-line per- 
turbation theory. Our method is applicable to any 
diagrammatic series expansion. We observe that in 
bare diagrammatic expansions, different possible contrac- 
tions for fermion operators typically sum up to deter- 
minants, substantially reducing the number of diagrams 
to be evaluated and ameliorating any minus sign prob- 



lem. In bold expansions the determinant structure is 
lost and the important question is whether the loss is 
offset by the physics gained from the partial resumma- 
tion. As a nontrivial example we apply the method 
to the one-impurity Anderson model, for which ana- 
lytical resummation techniques are well established and 
are believed to capture much of the physics, and a 
substantial body of numerical work exists for compar- 
ison. The analytical resummation techniques are the 
non-crossing approximatio n^^d^ ("NCA") and the one- 
crossing approximationi^ ( "OCA" ) , where the names re- 
fer to specific topological features of the diagrammatics 
which we discuss in more detail below. We formulate the 
bold-line expansion about these analytical resummations, 
delineate the regimes in which it is useful, and use the 
techniques to resolve a long-standing question concerning 
the form of the electronic spectral function near the edge 
of the Mott gap. As a by-product we determine the range 
of applicability of the NCA and OCA approximations for 
spin degeneracy N = 2. 




FIG. 1. (Color online) Typical diagrams arising in the hy- 
bridization expansion of the partition function of the An- 
derson impurity model: the four local states are described 
by wavy lines (|0), | t4-))- light, red (dark, blue) dashed lines 
(I t)> I i))- Light, red (dark, blue) solid lines denote hybridiza- 
tion functions, and empty (filled) circles local annihilation 
(creation) operators, a: NCA diagram, b: "bold" propagator 
(double dashed line) which resums NCA diagrams including 
a. c: diagram with crossing lines, d: Diagram in expansion 
which resums diagrams including c. In BoldOCA, all of dia- 
gram c is contained in the single bold line b. 

The rest of this paper is organized as follows: in SecHIl 
we outline the general features of our bold-line method. 
In Sec, mil we present the specific formulas needed for the 
application of the method to the single-impurity Ander- 
son model. In Sec. IIVI we present metrics demonstrating 
the successes and limitations of the method and in Sec.lVl 
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we show results for the electron Green's function, demon- 
strating in particular that the greatly improved accuracy 
of the method allows a definitive statement about the 
structure of the spectral function at the edge of the Mott- 
Hubbard gap. Section I VI I is a summary and discussion 
of future prospects. 



II. METHODS 

Bold-line expansions of many-body problems require a 
modification of the diagrammatic Monte Carlo technique 
(these modifications were not necessary in the single- 
particle problem studied in Refs. H and Conven- 
tional diagrammatic Monte Carlo is typically formulated 
in the configuration space consisting of all diagrams con- 
tributing to the partition function. We refer to this space 
as partition function space and denote it by Cz- In 
a conventional expansion the space Cz suffices because 
Grcen's-function diagrams are generated by breaking a 
propagator line in partition function diagrams. However, 
in a bold-line expansion the configuration space must be 
extended because there exist Green's-function diagrams 
which are not obtained by breaking lines in a bold ex- 
pansion of the partition function. We therefore employ 
a "Worm" algorithmi^ii^ Our calculation is formulated 
in an extended configuration space Cw consisting of the 
union of partition function [Cz) and Green's function 
(Cg) space: Cw = Cz U Cg, with a weight w(x) of a 
given configuration x given by wz if x S Cz and wq if 
x S Cg- In this extended space we define a new partition 
function 



= Z + ?; / G 



(1) 



where / denotes a sum over all the elements of the G 
matrix and the parameter 77 is in principle arbitrary, but 
in practice should be chosen so that the Z and G parts of 
W give comparable contributions; the calculations pre- 
sented here use 77 = 0.15. 

A measurement of a component ah of G, {Gah)W; in 
the extended space Cw is proportional to the physical 
Green's function but is wrongly normalized. The correct 
normalization is obtained after division by the partition 
function, also measured in the space Cw i i-e. by {?i'^)w = 
X]x3x ''^itli 3x = 1 if X g Cz and 3x = otherwise. Thus, 



ir< \ _ '\Gah)w 

K^ah) - -JTz 



w 



(2) 



The sum over all terms in Cw is performed using a 
diagrammatic Monte Carlo method: diagrams in Cw are 
generated, accepted, or rejected stochastically by insert- 
ing and removing local operators and hybridization lines 
according to their contribution to Zw and integrated 
stochastically, in analogy to Refi^. This summation is 
exact if all bold diagrams are included. 



While diagram-generating procedures are model de- 
pendent, they generically consist of the insertion or re- 
moval of operators and reconnection of parts of diagrams. 
The proposed insertion (removal) of operator tuples (typ- 
ically pairs) raises (lowers) the diagram order, changing 
a configuration x to y and is accepted with probability 



Tj^acc 



i?x 



min(l,i?xy) 
^;(y)WP^°P 



R 



yx J 



(3) 
(4) 



where W^y°^ denotes the proposal probability of an up- 
date. Transitions between configurations in Cg and Cz 
are performed by inserting or removing Green's function 
operators or propagator lincsinto (from) a partition func- 
tion configuration. Typical diagram generating proce- 
dures can produce terms which are already included in 
the bold resummation; one must test each generated di- 
agram to make sure that it is not already included, but 
the computational cost of the test is typically negligible. 



III. APPLICATION: SINGLE-IMPURITY 
ANDERSON MODEL 

As a stringent test of the bold-line resummation meth- 
ods wc consider one of the best-studied nontrivial models 
in condensed-matter physics, namely, the single-impurity 
Anderson model. This model represents the physics of a 
magnetic impurity in a metal and is also important as an 
auxiliary problem in the "dynamical mean- field" approx- 
imation to the properties of models of correlated electron 
materialsii For this model, widely used partial resumma- 
tion methods are available. In addition, numerically ex- 
act results obtained by the hybridization expansion (CT- 
HYB) quantum Monte Carlo method^ are available for 
comparison. 

The Anderson model describes a correlated site cou- 
pled to a bath of free electrons; the Hamiltonian is 



bath 



loc 



(5) 



with 



per per 

Hioc = ^ -~^J.ncr + U nfUi, 

a 

where the a^a label a continuum of "bath" operators with 
dispersion tp, and Ca arc local operators at energy /i with 
interactions [/, hybridizing with the bath with strength 
Vp. In the original formulation of the Anderson model, 
the parameters £d and C/ are properties of the magnetic 
impurity and Ep and Vp are properties of the host metal. 
In the dynamical mean-field context the impurity model 
is an auxiliary problem used to provide information about 
a lattice model of interacting electrons. In this case U is 
the on-site interaction of the lattice model and the £d, £p 
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and Vp are fixed by a self-consistency condition, as dis- 
cussed in Ref. 0- These details are not important for the 
formulation and application of the bold-line expansion. 
However, our specific results are obtained using parame- 
ters arising from the dynamical mean- field solution of the 
one-orbital Hubbard model at various particle densities 
and interaction strengths. 

We study the model as an expansion in V about the 
atomic limit in which the local states are decoupled from 
the bath. Alternatively, a resummation of interaction di- 
agrams around the free limit could be designed; we do 
not consider this here. The expansions may be repre- 
sented pictorially by time-ordered diagrams such as those 
shown in Fig. [1]. The propagation in atomic states 
is represented by different lines (wavy and dashed): in 
the eigenbasis |0), | t), I ti) of i?ioc the energies are 
ii^ioc — 0, —fi, — /i, U — 2/i, respectively, and the corre- 
sponding (bare) propagators are e""^^'"". The hybridiza- 
tion vertices are indicated by solid and open circles while 
the propagation of electrons in the bath is denoted by 
solid lines which represent the hybridization function 



Mr) = E ■ 



v*v 



T > 
T < 0. 



(6) 



to diagrams that contain no non-crossing parts (i.e., we 
use the CT-HYB method to generate diagrams but do 
not sample those [such as Fig.[IJc)] that have a part con- 
tained in NCA, and we replace the atomic propagators 
Gj'^-^ by NCA propagators G^^^^^). In order to sample 
around the first-order OCA, we limit ourselves to dia- 
grams without one-crossing subsegments [thereby reduc- 
ing Fig. [TJc) to Fig. [TJ6)]. Testing if a diagram needs to 
be sampled is a (cheap) 0(k) operation: in BoldNCA di- 
agrams with hybridization lines a^flp that span no other 
local operators need not be sampled. In BoldOCA, ad- 
ditionally, diagrams of the type aj^OgOpaJ are excluded. 

The CT-HYB diagram-generating procedure is ergodic 
(generates all diagrams); while we do not have a formal 
proof that it remains ergodic when restricted to bold dia- 
grams we have been unable to find a counterexample and 
have extensive numerical evidence that all bold diagrams 
are generated. 

An independent discussion of a bold expansion about 
the NCA limit was given in Ref. [l^. The extended config- 
uration space Cw was not employed. A bold expansion 
as defined here was not implemented but an evaluation 
of all diagrams up to 5th order, with time integrals eval- 
uated by a Monte Carlo method, was reported. 



Accurate and efficient numerical methods exist for eval- 
uating the direct hybridization expansioui^ 

The resummation methods we consider are the NCA 
(Refs. dH and ^ and the OCAJ^ The NCA resums 
all segments containing no crossing fermion lines [e.g.. 
Fig. [1] (a)] into a renormalization of the propagator of 
the atomic state [e.g.. Fig. [T](6)], which becomes G\j^ = 



equations 



S|j^] with Yj\j^ given by the self-consistent 



S|o)(r) = G|t)(r)At(r) + G,;) (T)A;(r), (7) 
S|.)(t) - G|o)(T)A,(-r) + G|t4)(T)A_.(T), (8) 
S|n)(T) = G|t)(r)A;(-r) + G,^) (t)A^(t). (9) 

First-order OCA self-energies arc given by equations in- 
volving additional crossing hybridization lines, and the 
full OCA equations are obtained by resumming vertex 
equations. Note that the projection techniques used in 
analytical NCA calculation a^^'^^ are not needed here. 
Multi-orbitalii and cluster models can be studied using 
the same operations with A, G, and S in matrix form. 
The hybridization expansion of the partition function 

is 



-''^-IVc,,(TV,)c],(Ti)...c,,(ri)c],(rO 



(10) 



X Tr, 



The hybridization expansion algorithm^iil (CT-HYB) 
enables a direct numerical sampling of this series. To 
sample a bold-line expansion around the NCA we restrict 



IV. RESULTS: METRICS 
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FIG. 2. Probability p for a diagram of Z to contain k spin-up 
hybridization lines, computed for Anderson impurity model 
with semicircular density of states and hybridization function 
fixed by dynamical mean-field self-consistency condition with 
parameters 77 = 0.15 n = 1, jSt = 10, at (7-values indicated. 
Full symbols: BoldNCA. Empty symbols: BoldOCA.. 

In this section we present some basic results which 
illustrate the strengths and weaknesses of the bold ex- 
pansion. Our results were obtained using the Ander- 
son model corresponding to the dynamical mean-field 
approximation to the single-orbital Hubbard model on 
a Bethe lattice.— This model is specified by an inter- 
action strength J7, a hopping parameter and a car- 
rier concentration n. At carrier concentration n ~ 1 
the dynamical mean-field approximation to this model 
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has a metal-insulator transition at a critical interaction 
strength U 5t (at the temperatures we study )i and we 
shall see that the bold expansion behaves very differently 
in the insulating and metallic phases. 

An important metric is the weight p{k) of terms at or- 
der k in the expansion for Z. This is shown in Fig. [2] 
for different interaction strengths U at inverse tempera- 
ture 13 = 10/t. The curves are characterized by a non- 
negligible weight at zero expansion order and a tail which 
decays approximately exponentially at high perturbation 
order. If the NCA (or first-order OCA) exactly solved 
the model p(fc = 0) would equal 1 and the tail would be 
absent. 

A clear difference between the curves is seen. For 
^ 5 p(fc = 0) is very close to unity and the weight 
in the tail is very small, while for U < 5t p{k = 0) be- 
comes smaller and the contribution of the non-NCA (non- 
OCA) diagrams becomes important. These differences 
arise from a difference in physics. At this temperature 
the model is in a gapped insulating phase for U > 4.5i 
and is in a gaplcss metallic phase for U < 4.5t. Clearly 
the NCA and first-order OCA are very good approxima- 
tions to the insulating phase and poor approximations 
to the metallic phase. A bold expansion around the free 
{U = 0) limit (not considered here) presumably behaves 
differently: while the weakly correlated metal would be 
captured accurately the local physics of the Mott insula- 
tor would be difficult to reach. 

We have also studied these histograms of p{k) as a 
function of temperature (not shown). The temperature 
dependence is negligible in the insulating phase. In the 
metallic phase, as the temperature is decreased, the value 
p{k = 0) decreases and the weight in the tail increases; 
as r — > p{k = 0) — > reflecting the failure of NCA and 
OCA to adequately describe the physics in the T —> 
limit of the metallic phase. Note that the weight of di- 
agrams as a function of expansion order as expressed by 
the long tail in Fig. [3] only decays very slowly and that 
basing the algorithm on a more sophisticated resumma- 
tion (here, first-order OCA instead of NCA) does not 
change the decay of the tail. 

Another important metric for an expansion of an in- 
teracting fermion problem is the average sign. Figure [3] 
shows that the bold expansion suffers from a sign prob- 
lem. The upper panel demonstrates that in the insulating 
regime (U > 4.5, at this T) the expansion around Bold- 
NCA shows convergence at very low order and configura- 
tions with negative sign (which occur at higher expansion 
order) give a negligible contribution. Thus, in this case 
the loss of the determinant structure is compensated by 
the much better starting solution. However, in the metal- 
lic case the starting solution is less good and a severe sign 
problem arises. The steep drop around U/t = 4 marks 
the departure from the insulator and the failure of NCA 
and OCA. The lower panel shows a similar behavior of 
the sign as a function of doping at a strong interaction. 
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FIG. 3. Upper panel: expectation (over Cw) of sign as a 
function of interaction U , at carrier concentration n = 1 per 
site and temperatures indicated. Lower panel: expectation 
value of sign as a function of carrier concentration 1 — a; at 
U = 6t and temperatures indicated. 



V. RESULTS: ELECTRON GREEN'S 
FUNCTION 

Figures |4] and [5] compare the imaginary time Green's 
function obtained by BoldNCA (BoldOCA results are in- 
distinguishable on the scale used in Figs. |3]and[S]) simu- 
lations to the CT-HYB results and analytical NCA and 
OCA results. Comparable computational resources are 
invested in the bold and CT-HYB calculations. The large 
differences between the NCA and OCA and the Monte 
Carlo curves in the upper panel of Fig. Ushow that for the 
weakly correlated metallic phase the initial starting point 
is poor and the bold expansion is not useful in practice. 
The lower panel shows that in the moderately correlated 
"bad metal" case the starting point is closer to the ex- 
act answer and the bold and CT-HYB results are very 
close at temperature T = t/10. The rapid decrease of 
sign with T (cf. Fig. ^ shows that for U/t = 4, /3i w 15 
the bold method will become significantly less efficient. 
The upper panel of Fig. [5] shows that the same situa- 
tion is obtained in the strongly correlated, lightly doped 
case. Finally, the lower panel of Fig. [S] demonstrates the 
clear superiority of the bold methods in the insulating 
case, where the noise in the center of the imaginary-time 
interval is very substantially reduced. 

In standard CT-HYB simulations, spectra of insulat- 
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FIG. 4. Green's function of Anderson impurity model with 
semicircular density of states and hybridization function fixed 
by a DMFT self-consistency condition, calculated using NCA, 
OCA, CT-HYB and BoldNCA (bold OCA would be indistin- 
guishable), I3t = 10, starting from converged and accurate^fl 
CT-HYB results, at interactions and dopings indicated. Up- 
per panel: Fermi liquid case U/t = 1. Lower panel: Corre- 
lated Metal U/t = 4. 



FIG. 5. Green's function of Anderson impurity model with 
semicircular density of states and hybridization function fixed 
by a DMFT self-consistency condition, calculated using NCA, 
OCA, CT-HYB and BoldNCA (bold OCA would be indistin- 
guishable), pt = 10, starting from converged and accurate^S 
CT-HYB results, at interactions and dopings indicated. Up- 
per panel: Doped Mott insulator at filling n = 0.95 and in- 
teraction U/t — 6. Lower panel: Mott insulator, U/t = 6, 
n = 1. 



ing systems are hard to obtain because of large (relative) 
errors in the mid-range of the imaginary-time interval 
and so the questions of the precise value of the insulating 
gap and the form of the above-gap structure have been 
discussed extensively in the literature; for recent work 
and references, see Refs. [2TI and [2^ Figure [5] shows spec- 
tral functions obtained by maximum entropy analytical 
continuation of the bold-line, CT-HYB, and analytical 
NCA and OCA approximations. The relative errors in 
the basic CT-HYB data lead to a substantial smearing 
of the gap edge features. Also shown is the continuation 
of high-precision CT-HYB data obtained by expending 
64 times more computer resources. The bold expansions 
and the high-precision CT-HYB data essentially agree on 
the gap value and the form of the spectral function near 
the gap edge. (The differences at higher frequency reflect 
the intrinsic sensitivity of analytical continuation to very 
small differences in data.) 

Our results firmly establish that the spectral function 
in the paramagnetic insulating phase is characterized by 
a sharp peak at the gap edge, and wc obtain more precise 
values for the insulating gap. The accuracy is confirmed 
by the lower panel of Fig. [B] which presents the differ- 



ence AG between the measured G(t) and G(t) back- 
continued from A{uj). The analytical approximations 
show clear deviations, in particular a larger gap. The 
difference in gap value can be seen directly as a difference 
in the imaginary-time data in the lower panel of Fig. [S] 
G(t) calculated using the two analytical methods falls be- 
low the numerically exact results in the imaginary-time 
range 2 < rt < 5. 



VI. CONCLUSIONS 

In conclusion, we have presented a "bold" diagram 
method and applied it to an expansion around the non- 
crossing and one crossing approximations to the Ander- 
son impurity model. We have also applied the algorithm 
to the Kondo limit of the Anderson model without a self- 
consistency condition, with results (not shown) very sim- 
ilar to those described here. The NCA and OCA approxi- 
mations are believed to become exact in the limit of large 
iVi^ but their accuracy for physically relevant N has 
been established here. BoldNCA is general, numerically 
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— CT-HYB 

— Continued NCA 

— Continued (full) OCA 

— BoldNCA 

— BoldOCA 

-- CT-HYB, X 64 




FIG. 6. Upper panel: Analytic continuation of Green's func- 
tions obtained from different solvers from DMFT of half-filled 
Hubbard model at f/ = 6t and /3t = 7. Lower panel: difference 
AG(r) between G(r) obtained from CT-HYB, BoldNCA, and 
BoldOCA measurements and G(t) reconstructed from the an- 
alytically continued A{u}). 



exact, and easily extensible to multiple orbitals or cluster 
calculations, where a severe sign problem appears also in 
the hybridization expansion and the effect of reducing the 
expansion order will be most pronounced. We expect the 
main region of applicability to be in this well as 

to nonequilibrium problems, where the oscillating phase 
severely limits the applicability of non-bold methodsi^ 
Extensions to non-equilibrium and multi-orbital calcula- 
tions and to the resummation of vertices are currently 
under way. Diagrammatic random phase approximation 
('RPA') resummation is also essential for the study of 
screening and polarization; exploration of the methods 
discussed here may be fruitful in this context. 



ACKNOWLEDGMENTS 

We thank X. Wang for providing continuations shown 
in Fig. m E. G. and A. J. M. are supported by NSF 
under Grant No. DMR-0705847, and acknowledge par- 
tial support from NSF Grant No. PHY05-51164. Our 
programs are based on the ALPS library.— Calculations 
were performed on the Brutus cluster at ETH Zurich. 



^ N. V. Prokof'ev, B. V. Svist unov, 
JETP Letters 64, 911 (1996)1 



and I. S. Tupitsyn, 



N. V. Prokof'ev and B. 
Phys. Rev. Lett. 81, 2514 (1998) [ 
A. N. Rubtsov, V. V. Savkili^ and A 
[Pi^As. Rev. B 7 2, 035122 (2005) 
* P. Werner, A. Comanac, L. de' 
|Phys. Rev. Lett. 97 076405 (2006) 



V. 



Svistunov, 



I. Lichtenstein, 



Medici, et al, 



E. Gull, P. Werner, O. Parcollet, and M. Troyer, Europhys. 
Lett. 82, 57003 (2008). 

P. Werner , T. Oka, and A. J. MiUis, 

Phys! Rev. B 79, 035320 (2009) , 
A. Georges, G. Kotliar, W. Krauth 



and M. J. Rozenberg, 



Rev. Mod. Phys. 68, 13 (1996)| 



L. Miihlbacher and 

Phys. Rev. Lett. 100, 176403 (2008)| 



E. 



Rabani, 



N. Prokof'ev and B. Svistunov, "Bold diagrammatic monte 
carlo: A generic technique for polaron (and many-body?) 
problems," (2008). 

N. Prokof'ev and B. 



Phys. Rev. Lett. 99, 250201 (2007) 



H. Keiter and J. 

|J. Appl. Phys. 42, 1460 (1971)| 



C. 



Svistunov, 



Kimball, 



N. Grewe and H. Keiter, |Phys. Rev. B 24, 4420 (1981M 
T. Pruschke and N. Grewe, |Z. Phys. B 74, 439 (1989) 



E. Burovski, N. Prokof'ev, B. Svistunov, and M. Troyer, 
New J. Phys. 8, 153 (2006). _ 
N. E. Bickers, Rev. Mod. Phys. 59, 845 (1987) 
A. C. Hewson, The Kondo Problem to Heavy Fermions 
(Cambridge University Press, 1993). 

P. Werner and A. J. Millis, 



Phys. Rev. B 74, 155107 (2006) 



K. Haule, C.-H. Yee, and K. Kim, "Dynamical mean- 
field theory within the full-potential methods: Electronic 
structure of ce-115 materials," (20091 ■ larXivV09 07.0195vli 
M. Imada, A. Fujimori, and Y. Tokura, 



|Rev. Mod. Phys. 70, 1039 (1 998) 
2° E. Gull, P. Werner, "ST Mi llis 
|Phys. Rev. B 76, 235123 (2007)1 

21 



and M. Troyer, 



X. Wang, E. Gull, L. de' Medici, M. Capone, and A. J. 
Millis, Phys. Rev. B 80, 045101 (2009) 
I. S. Krivenko and A. N. Rubtsov, "Analytical continuation 
of quantum monte carlo data: Optimal stochastic regular- 



ization approach," (2006), arXiv:cond-mat/0612233 



A. Albuquerque, F. Alet, P. Corboz, et i 
Mater. 310, 1187 (2007). 



:., J. Magn. Magn. 



